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We introduce the spectrum bifurcation renormalization group (SBRG) as an improved application 
of the excited-state real space renormalization group (RSRG-X) for a class of qubit models. Starting 
from a disordered many-body Hamiltonian in the full many-body localized (MBL) phase, SBRG 
flows to the MBL fixed-point Hamiltonian, and generates the local integrals of motion and the matrix 
product state representations for all eigenstates. The method is applicable to many interacting spin 
and fermion models in the full MBL phase. As a Hilbert-space preserving RG, SBRG also generates 
an entanglement holographic mapping, which duals the MBL state to a fragmented holographic 
space. 
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The phenomenon of many-body localization ^BLH"^ 
has attracted much research interest recentlyPHSSl R is 
the generalization of the Anderson localizatiorP^to inter¬ 
acting syst ems, and can be described as the Eock space 
localizatioiP^^^of finite energy density states in the pres¬ 
ence of strong disorder. An MBL system cannot serve as 
its own heat bath, and thus violates the eigenstate ther- 
malization hypothesis (ETH) Erom man y aspe cts, 
the MBL excited states are like ground states.lI^IHl Por 
example, the entanglement entropy (EEp3 follows the 
area-law scalin^^ in the MBL state,li^*^ in contrast to 
the volume-law scalin^^SHSl the ETH state. In the full 
MBL system,!^ all energy eigenstates are localized and 
have area-law entanglement entropy, which implies the 
existence of matrix product state (MPS^PtM! 

represen¬ 
tations for all eigenstates. 

A natural question is whether there is an efficient 
way to find (or approximately find) the MPS represen¬ 
tation for every eigenstate in the full many-body spec¬ 
trum of a MBL system? In fact, there have been 
several nice approaches trying to answer this ques¬ 
tion, either by matching the unitary circuit of exact 
diagonalization,!^ or b y the spectral tensor network for 
local integrals of motion,l2Sl^or by variational tensor net¬ 
work optimization.!^ Here we would like to tackle the 
problem from a different angle using the renormalization 
group (RG) approach. For exa mple, the density matrix 
renormalization group (DMRGr^^^S been shown to 
be a powerful and successful method to find the MPS 
representations for ground states. For full MBL systems, 
DMRG can be generalized to target highly-excited states 
in the many-body spectrum as well, known as DMRG- 
^mHUl DMRG-X is a non-perturbative high-accuracy 
approach to find the MPS representation for an individ¬ 
ual excited state at a given energy-density. However it 
is still quite expansive to obtain all the eigenstates using 
DMRG-X. We would like to propose a less accurate but 
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more efficient RG approach to resolve all eigenstates in 
the MBL spectrum simultaneously. Obviously, we must 
take advantage of the strong disorder nature of the MBL 
system in order to target the full spectrum. The real- 
space renormalization group (RSRG^PEH 

(also known 

as the strong disorder renormalization group) is such an 
RG scheme that is designed to work with strong disor¬ 
der. The original proposal of RSRG is ground state tar¬ 
geting, and it was rec ently extended to target all excited 
states as RSRG-X.l^^lEoMl RSRG-X finds all eigenstates 
by transversing the spectrum branching tree explicitly. 
However for a specific class of qubit models (such as quan¬ 
tum Ising models), where the spectrum branches at each 
RG step to two subspaces with identical dimensions, it is 
then possible to implicitly encode the spectrum branch¬ 
ing by an emergent qubit degree of freedom, and formu¬ 
late the RG as a flow of the entire many-body Hamilto¬ 
nian without branching to any specific energy subspaces. 

The idea of formulating RSRG-X as a Hamiltonian flow 
has appeared in several recent works Ref. lTBHTSl In this 
work, we will introduce a systematic RG scheme along 
this line of thought, as outlined in the following. First 
we pick out the leading energy scale (highest frequency, 
largest coupling) term Hq in the Hamiltonian iJ, and 
rotate the Hamiltonian H to the diagonal basis of Hq 
by a local unitary transformation. Next we eliminate 
the terms that anti-commute with the leading term Hq 
by second order perturbation, and obtain a new Hamil¬ 
tonian. Then we turn to the next leading energy scale 
term in the new Hamiltonian and repeat the above RG 
procedure, until all degrees of freedom in the system are 
accounted for. Finally, we collect the unitary transforma¬ 
tions that have been performed along the RG flow, and 
arrange them into a tensor network, then the MRS rep¬ 
resentations for all eigenstates are found and encoded in 
the tensor network. Because each RG step bifurcates the 
spectrum until the Hamiltonian is eventually diagonal¬ 
ized, we decided to name the RG scheme as the spectrum 
bifurcation renormalization group (SBRG). 

SBRG is a Hilbert space-preserving RG scheme, keep¬ 
ing track of the flow of the whole many-body Hamilto¬ 
nian without truncating the Hilbert space to any specihc 
branch of the spectrum. It can be considered as an im¬ 
provement of RSRG-X, in the sense that SBRG allows 
generation of new terms under the Hamiltonian flow, and 
also the history dependence of the spectrum branching 
is not explicitly realized but implicitly encoded in the 
ultraviolet-infrared (UV-IR) mixing terms in the Hamil¬ 
tonian. Therefore SBRG can be applied to strong inter¬ 
acting models beyond the limitation of a closed-form RG. 
Like RSRG-X, SBRG also targets the full many-body 
spectrum rather than just the ground state, and focuses 
on small frequency (small energy difference) rather than 
the low absolute energy. However, SBRG also suffers 
from the same limitation as RSRG that it is a basis de¬ 
pendent approach, which means local unit ary rotations of 
the Pauli basis will affect the RG flow.E^ This limitation 
may be overcome by resolving the local basis rotation be¬ 


fore identifying the leading energy scale, however we will 
leave the possible improvements for future study. 

Another motivation comes from the recent research 
effort to interpret the RG transformation as a realiza¬ 
tion of the holographic duality,^SlH2ni examples include the 
multiscale entanglement renormalization ansatz (MERA) 
networtPS^^, cMERA^^^lt^, ab initio holographjPSl ect. 
In particular, Ref.lHUHH proposed the entanglemenipi^l 
holographic mapping (EHM) as a unitary mapping be¬ 
tween boundary (physical) degrees of freedom and bulk 
(emergent) degrees of freedom for every Hilbert space¬ 
preserving RG. Here we apply the same idea to SBRG 
for the MBL system, and it turns out that the emer¬ 
gent degrees of freedom are just the emergent conserved 
quantitied^mss] MBL fixed-point Hamiltonian. 

The emergent conserved quantities are identified by 
SBRG progressively as the leading energy scale. The 
spectrum bifurcation at each RG step is controlled by 
the emergent conserved quantities (like controlled-gates 
in the quantum circuit). In the end of the RG flow, the 
MBL fixed-point Hamiltonian will emerge. By collecting 
the unitary transformations that have been performed 
to the Hamiltonian during the RG flow as matrix prod¬ 
uct operators (MPO),l^one can reconstruct the approx¬ 
imate tensor network representation of the full many- 
body spectrum, which encodes the approximate MPS 
representation of each MBL eigenstate. An important 
observation of this work is that the EHM circui t of SBRG 
can be approximated by a Glifford circuiiP^^^ which en¬ 
ables efficient calculation of many physical properties of 
the MBL system. 

The paper is organized as follows. In Section we 
start by introducing SBRG algorithm under generic set¬ 
tings, and discuss the limitations and strengths of the 
method. In Section in we apply SBRG to the disor¬ 
dered quantum Ising model, and benchmark the energy 
spectrum and eigenstates obtained by SBRG. We will 
demonstrate that SBRG flows towards the strong disor¬ 
der limit in the MBL phase. The RG flow also generates 
a Glifford circuit which encodes the (approximate) MPS 
representations of all eigenstates. In Section |IV[ we in¬ 
terpret the Glifford circuit as an EHM tensor network. 
We investigate the properties of the emergent conserved 
quantities and conhrm the area-law entanglement in the 
MBL phase. Finally we will discuss the locality in the 
holographic bulk, and the collapse of SBRG near ther- 
malization. 


II. SPECTRUM BIFURCATION 
RENORMALIZATION GROUP 

A. SBRG Algorithm 

The basic idea of SBRG is to progressively identify 
emergent conserved quantities by block-diagonalization 
of the leading energy scale, and eliminate the block-off- 
diagonal terms in the Hamiltonian by second order per- 
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turbation. The most general setting is to start from the 
following qubit model Hamiltonian 

^ = = (1) 
[m] [a«] 

where = 0,1,2,3) denotes the Pauli matrix acting 

on the zth qubit. = cr/"iM 2 M 3 - = cr^i (gjcr'^a 0 • • • 
is a short-handed notation for the direct product of Pauli 
matrices, or called Pauli operators. The Hamiltonian is 
simply a sum of the Pauli operators with real coeffi¬ 
cients over all Pauli indices [/i], which are randomly 
drawn from independent distributions. Here we have as¬ 
sumed the Hilbert space dimension is a power of two as 
2^, such that any Hamiltonian in the Hilbert space can 
be written as a linear combination of the Pauli opera¬ 
tors crl^l. The qubit model is quite general. It allows 
multi-qubit non-local interactions on a generic lattice, 
and can describe a large class of spin and fermion sys¬ 
tems in all dimensions (fermion models can be mapp ed 
to qubit models by Jordan-Wigner transformatiorpdil). 

Starting from the qubit model Eq. Q, SBRG algo¬ 
rithm goes as follows. First, we pick out the leading 
energy scale term in the Hamiltonian, which amounts to 
selecting the term with the maximal coefficient 

I (in its absolute value) among all terms in the Hamil¬ 
tonian H, and denote it as 

iJo = (2) 

where {hsl = = maxj^j represents the lead¬ 

ing coefficient (the meaning of the subscript 3 will be 
evident later). Because each Pauli operator has nor¬ 
malized eigenvalues ± 1 , the energy scale of each term 
is only determined by the coefficient in the front. 
At this point, we need to assume that there is only one 
unique term with the leading energy scale, and all the 
other terms have energy scales sufficiently less than the 
leading one. This assumption can be justihed in the 
strong disorder limit, but will break down for uniform 
(translational invariant) systems. 

Then we block-diagonalize the leading term by 

a Clifford group rotation R such that 

^ (3) 

and at the same time transform all the other terms in 
the Hamiltonian by the same rotation, such that H —>■ 
B)HR. Here denotes cr^ times the rest of the 

identity matrices. As a technical note, although we write 
the diagonal form as for theoretical formulation 

of the RG scheme, in practice the operator does not 
need to be swap to the first qubit literally. Its action can 
remain in the local support of the original Pauli operator 
^Mmax^ such that i? is a local unitary transformation (see 
Fig.[3llater for an explicit illustration, and also Appendix 
for implementation details). The key observation is 
that R is not a generic unitary transformation, but an 
element in the Clifford group which rotates one Pauli op¬ 
erator to another. So R can be easily found for any given 


(see Appendix A1 for the algorithm), and can be 


applied to the Hamiltonian efficiently. This is related to 
the Gottesman-Knill theorenP^ that Clifford circuits can 
be simulated efficiently on a classical computer. 

As we block-diagonalize the leading energy scale, Hq 
becomes Hq = and the many-body spectrum 

bifurcates to the high-energy A ~ I/ 13 I and the low-energy 
E ~ — 1 / 13 1 sectors (blocks). We must reduce other terms 
in the Hamiltonian into either the higher or the lower 
energy sectors. Thus we classify the terms by their com¬ 
mutativity with Hq as 


H = Ho + A + E, (4) 


where A are terms that commute with Hq, i.e. AHq = 
HqA', and E are terms that anti-commute with Hq, i.e. 
EHq = —HqYj. All terms must fall into these two classes, 
because Hq contains only one single term of a Pauli oper¬ 
ator (not the sum of several terms). Given Hq = 
in the i?-rotated basis, one can see that A rests in the 
diagonal block as a combination of l and , 

and E rests in the off-diagonal block as a combination of 
and The diagonal terms A is left untouched, 

and are passed down with the Hamiltonian to the next 
step of SBRG. The off-diagonal terms E must be renor¬ 
malized by 2 nd order perturbation, which corresponds 
to the unitary transformation H —>■ HS (known as the 
Schrieffer-Wolff transformatiorP^ with 




(5) 


A2 


for 


carried out to the 2nd order of h'^ (see Appendix 
derivation). The resulting effective Hamiltonian within 
the diagonal blocks reads 


H = Hq + A- Eht^E 


= Hq 


■A 


( 6 ) 


2hl 


HqE^ 


Under the Schrieffer-Wolff transformation S, the Hamil¬ 
tonian is block-diagonalized (to the 2 nd order), as one 
can check that Hq'E'^ now commutes with Hq. Since new 
terms are generated under the 2 nd order perturbation, 
the number of terms in the Hamiltonian will presumably 
grow in this step. In practice, small terms can be trun¬ 
cated to control the growth rate (see Appendix A 2 for 
the truncation scheme). 

So finally, the new Hamiltonian takes the form of 




where the first qubit is acted by with A = 0 or 3 only, 
and the remaining qubits are acted by the generic with 
^ = 0,1,2, 3. The leading term h 3 cr^l° ' 1 is singled outP^ 
and is then ascribed to the effective Hamiltonian (which 
will not be touched in the later RG steps). The oper¬ 
ator is also identified as an emergent eonserved 
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quantity, because its energy scale is so high (compared 
to its local neighbors) that it is unlikely to flip as we 
zoom in the spectrum towards the low frequency limit. 
The remaining terms i7res = Sa[/j] have the 

same form as the qubit model Eq. 0 that we started 
with. They will enter the next SBRG step. In the next 
step, we continue to pick out the leading energy scale 
term in i?resj diagonalize it by a Clifford group rota¬ 
tion to the second qubit as and reduce the off- 

diagonal terms by the 2nd order perturbation. Because 
i/res commutes with it can be diagonalized by the 

(j3[0 "]_preserving unitary transformations, meaning that 
the previously identihed conserved quantity will 

not be affected by the later RG transformations. 

As we keep running SBRG procedure, the Hamiltonian 
will flow to the following generic form, which splits into 
two parts (where A = 0, 3 and /r = 0,1, 2, 3) 

H = HeS + i?res, 

[A] (8) 

[A].M 


spectrum by Monte-Carlo samplin^^ or thermally aver¬ 
ages over all branched^. In SBRG, we keep the spec¬ 
trum bifurcation structure implicitly in the Hamiltonian 
on the operator level, and keep track of the flow of the 
whole Hamiltonian in the full Hilbert space. 



FIG. 1: Flow of the many-body Hamiltonian under SBRG. 
Each RG step consists of a Glifford group rotation Rk followed 
by a Schrieffer-Wolff transformation Sk, which block diagonal¬ 
izes the Hamiltonian and identifies an emergent qubit. In the 
end, the Hamiltonian will flow to the diagonal form Heg. 


SBRG ends after the last qubit in the physical Hilbert 
space has been renormalized to the emergent Hilbert 
space. Then we are left with the effective Hamiltonian in 
terms of the emergent conserved quantities only 

HeS^ ( 10 ) 

[A=0.3] 


The effective Hamiltonian Hgg contains all the terms 
that have been fully diagonalized in previous RG steps, 
and the residual Hamiltonian i?res contains the remain¬ 
ing terms to be diagonalized in future RG steps. The 
Hilbert space is also naturally partitioned into the emer¬ 
gent Hilbert space acted by crt'*'! ((A = 0,3), and the 
physical Hilbert space acted by (^ = 0,1,2, 3). In 
each RG step, the emergent Hilbert space grows by one 
qubit, while the physical Hilbert space shrinks by one 
qubit. Gorrespondingly the terms in H,-es are progres¬ 
sively renormalized to H^g. Each emergent qubit is an 
emergent conserved quantity identified at a particular en¬ 
ergy scale, which controls the branching of the spectrum 
of that energy scale. Since the Hilbert space is not trun¬ 
cated, the information of spectrum bifurcation is kept 
in the [A] dependence of the = crt'''] 0 terms 

in Hj-es- If we specify an Ising configuration |t) for the 
emergent qubits, (Tjcrl'’'! |r) will acquire definite expecta¬ 
tion values ±1, and the coefficient in front of the 
remaining term will be determined for this specific 
branch of the spectrum given by the Ising configuration 


(r|i?res|T) = 

(9) 

hn] 

[A] 

In this way, a specific branching choice can be made. 
However, we will not make such an explicit branching 
choice in SBRG flow. This is in contrast to RSRG-X 
approach, in which one either visits each branch of the 


which is fully diagonalized. So SBRG can be consid¬ 
ered as an approximate approach to diagonalize a ran¬ 
dom many-body Hamiltonian in the strong disorder limit, 
as illustrated in Fig.[l] The MBL fixed point Hamilto¬ 
nian proposed in Ref IzTVE!?! also takes the same form as 
Eq. (10), and the key point there was that the emergent 
conserved quantities are localized in terms of the original 
physical qubits, which can be verified in our numerics 
later. 


B. Limitations and Strengths of SBRG 

SBRG is a generalized and improved application of 
RSRG-X for a class of qubit models whose spectrum bi¬ 
furcates at each RG step. We should admit that SBRG is 
not a generic solver of disordered many-body Hamiltoni¬ 
ans. The method has several limitations: (i) SBRG only 
works for Ising/Majorana-like qubit models which have 
the bifurcating spectral structure, (ii) Within the scope 
of qubit models, SBRG still suffers from the problem of 
Pauli-basis dependency, (iii) Even under a nice choice 
of the basis, SBRG only works in the strong disorder 
regime, and breaks down near the thermalization tran¬ 
sition. Limitations (ii) and (iii) are shared hy RSRG-X, 
where (ii) can be overcome by DMRG-XP^"^ and other 
MPO-based RGP^. We will explain these limitations in 
details as follows. 

The current scheme of SBRG is designed to work with 
the qubit model Eq. ([^ , and crucially relies on the as¬ 
sumption that there is only one non-degenerate leading 
energy scale in each RG step. It will not work, for ex¬ 
ample, with the models built from three-state quantum 
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rotors, or with SU(2) symmetric spin models which have 
more than one terms of the same leading energy scale and 
singlet-triplet spectrum branching (which is not bifurcat¬ 
ing). To run SBRG, we must assume the many-body 
spectrum to follow certain bifurcation structures in a 2^ 
dimensional Hilbert space, such that the RG can follow 
the bifurcation tree and discover the full spectrum pro¬ 
gressively. That is why this RG scheme is called “spec¬ 
trum bifurcation”. Of course, further improvements can 
be made to generalize the current RG scheme to work 
with non-bifurcating spectrums, however we will leave 
this direction for future research. 

SBRG is not a basis independent RG scheme. It re¬ 
quires an fine-tuned Hamiltonian where the leading en¬ 
ergy scale term at each RG step must be given by a sin¬ 
gle Pauli operator. A local change of basis could spoil 
this requirement. For example if turns out to be 
a leading energy scale term. Under an arbitrary two- 
qubit basis rotation, say it could be transformed 

to h cos + /i sin 0(7°^, then each of the terms is not 
the true leading energy scale anymore. Thus local ba¬ 
sis transformation could affect SBRG flow, although the 
physics should not be affected. So the Hamiltonian must 
be written in a nice basis such that the leading energy 
scale term is always represented as a single Pauli oper¬ 
ator (or at least approximately). In general, local basis 
rotations will introduce complicated correlations among 
the coefficients /i[^] in the qubit model Eq. Q. So if we 
assume the coefficients are all independent (uncor¬ 
related), then in the strong disorder limit, nearby Pauli 
operators will have very different energy scales. In this 
case, the strongest Pauli operator will be a good approx¬ 
imation of the true leading energy scale term locally. 

As a kind of strong disorder RG, SBRG only works 
in the MBL phase and breaks down near the MBL-ETH 
transition (see Fig.[^. Two things could happen when 
we approach the thermalization transition from the MBL 
side, (i) Uncontrolled number of terms could be gen¬ 
erated by the 2nd order perturbation in Eq. ®> , such 
that the terms in the Hamiltonian could grow exponen¬ 
tially (or even faster), which would crash SBRG program, 
(ii) The off-diagonal terms E become close in energy scale 
to the leading term Hq, such that the perturbative treat¬ 
ment is no longer valid. However if we keep away from the 
ETH phase, SBRG should work well in the MBL phase 
and on the marginal MBL boundarjl^ (also known as the 
quantum critical glas^^ between two MBL phases (see 
Fig.§. Although delocalization happens at the marginal 
MBL criticality (i.e. the MBL-MBL transitiorfi^IlSl)^ yet 
the growth rate of the Hamiltonian terms is still con¬ 
trolled by the MBL phases from both sides. Thus SBRG 
is still applicable to the marginal MBL system. 

Despite of the above limitations, SBRG inherited the 
high efficiency of RSRG. It is fast compared to DMRG- 
X, and can handle relatively large system size (typical 
system sizes used in this work is of 128 to 512 qubits). 
Of course, the efficiency is gained at the cost of losing 
the accuracy. However if our goal is to study the over- 


cn 

c 
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« 

a> 

■D 

o 

!2 
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FIG. 2: Schematic phase diagrams containing two different 
kinds of transitions: the MBL-ETH transition (in red) and 
the MBL-MBL transition (marginal MBL, in black). At weak 
disorder, eigenstates in a many-body system typically ther- 
malized. However as the disorder gets strong enough, they 
can be many-body localized. The MBL phase can be further 
divided into different quantum phases by various quantum 
orders. Different MBL phases are separated by the marginal 
MBL critical line. Whether the marginal MBL to ETH tran¬ 
sition happens at a finite disorder strength is still an open 
question. SBRG fails in and near the ETH phase, but works 
well at the MBL-MBL transition and in the MBL phase. 


all structure of the MBL spectrum other than individual 
states, SBRG turns out to be an efficient method. To im¬ 
prove the accuracy, it is possible to pass the MRS states 
generated by SBRG as the initial states to DMRG-X for 
further improvement. 

The major improvement of SBRG compared to the 
generic RSRG-X is to encode the spectrum branching 
history implicitly in the Hamiltonian flow. This allows 
SBRG to obtain the MBL fixed-point Hamiltonian in one 
run. Many physical properties of the fixed-point Hamil¬ 
tonian can be therefore investigated, including the scaling 
and distribution of the energy coefficients and the emer¬ 
gent conserved quantities (to be elaborated in the follow¬ 
ing sections). SBRG can also generate an unitary MPO 
that approximately diagonalizes the MBL Hamiltonian. 
It turns out that the MPO can be approximately express 
as a Glifford circuit which has a high computational effi¬ 
ciency. With the efficient MPO, a controlled holographic 
mapping of the entire many-body Hilbert space becomes 
possible, which provides us some geometric interpreta¬ 
tions of the entanglement structures in the MBL states. 
Finally as a numerical approach, SBRG allows new terms 
to be generated with the RG flow, which overcomes the 
limitation of closed-form RG (as RSRG-X was originally 
proposed), and enables us to study strongly interacting 
and higher dimensional MBL systems. 


III. APPLICATION TO QUANTUM ISING 
MODEL 

A. Beyond the Closed-Form RG 

To illustrate and to benchmark the SBRG scheme, we 
will take the disordered quantum Ising model (transverse 


marginal MBL 

MBL MBL 

(phase A) ,, (phase B) 

2X2 

ETH 


tuning parameter 
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field Ising model plus interaction) as an example. The 
model is defined on a ID spin chain, with both XX and 
ZZ coupling and the external field Z, 

H = -Y^ Jicrlcrl^i + + h^af. ( 11 ) 


The model can also be interpreted as the interacting Ma- 
jorana fermion chain under the Jordan-Wigner transfor¬ 
mation, 



( 12 ) 


The coefficients Ji, Ki and hi are random variables in¬ 
dependently drawn from uncorrelated distributions. In 
the following, we will switch between the spin and the 
fermion interpretations whenever which one is more con¬ 
venient. 
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FIG. 3: Illustration of a single RG step, with the leading 
energy scale term Hq in red, the block-diagonal terms A in 
green, and the block-off-diagonal terms S in blue. The Hamil¬ 
tonian is first transformed to the diagonal basis of Hq by a 
Clifford rotation R. At this step, one qubit (marked in red) 
is identified as the emergent conserve quantity. Then the 
block-off-diagonal terms Ji and Js are treated as perturba¬ 
tions, which generate effective couplings among the rest of the 
qubits via a Schrieffer-Wolff transformation S. 


Fig.|^. The leading energy scale term is identified 

as a new emergent conserved quantity, and is ascribed 
to the effective Hamiltonian. The spectrum is also bifur¬ 
cated to high- and low-energy subspaces, depending on 
whether the emergent qubit operator takes -1-1 or —1 
eigenvalues. The remaining terms are passed down to the 
next RG step. The difference with the non-interacting 
case is that the residual Hamiltonian in the physical 
Hilbert space no longer keeps the form of Eq. 0 , new 
terms are generated: including the interaction among 
more qubits like {JiJ 3 /K 2 )cr^°^^ and the Zeeman field 
in other directions like J 2 cr°^^°. Some of the terms may 
contain the operator tr® on the emergent qubit, which 
encodes the dependence of the spectrum branching. So 
in general, it is no longer possible to preserve the RG 
transformation in a closed form in terms of a few deci¬ 
mation rules. Therefore SBRG is in principle a numerical 
method, which can be used to explore the full MBL phase 
diagram of Eq. (11) in the strong interaction regime. 


However before running the numerics, let us mention 
several limits of the disordered quantum Ising model in 
Eq. (II) which are well understood, (i) In the Ki = 0 
limit, the model simply goes back to the transverse field 
Ising model, and corresponds to a free Majorana chain. 
When In J > In h, the model is in a spin glass (SG) phase 
of the Ising chain (or as a topologically nontrivial phase 
of the Majorana chain). When In J < In/i, the model is 
in a paramagnetic (PM) phase of the Ising chain (or as a 
topologically trivial phase of the Majorana chain). Both 
of them are MBL phases, and they are separated by a 
marginal MBL critical point at In J = In h, which is also 
known as the infinite randomness critical point, (ii) In 
the hi = 0 limit, one can make a basis rotation to rewrite 


Eq. (11) as 


H — 'y ( Ji^i O'i+l + Ki^i 

i 


(14) 


In the original proposaP^ of RSRG-X, the bond term 
Jiolalj^i and the site term hiaf are considered as leading 
energy scales, and the “interaction” term Kiaiai_^i is 
treated as perturbation. Now in SBRG, all terms are 
allowed to be the leading energy scale. Eor example, 
suppose K 2 is the leading energy scale in the following 
4-site problem, one step of SBRG goes as: 

H = -K2a°^^° - 

_ - ... 

A _ j^^022i 

_ . ( 13 ) 

A + 0300 'll >^3 1011 

which consists of a Glifford rotation R followed by 
a Schrieffer-Wolff transformation S (as illustrated in 


which can then be converted to two independent copies of 
free Majorana chains via the Jordan-Wigner transforma¬ 
tion as illustrated in Eig.|^ Therefore although Ki ^ 0, 
the model Eq. (14) is still a free fermion model, and have 
a marginal MBL criticality at In J = In AT, which can be 
considered as the doubled version of the In J = In h crit¬ 
icality. (iii) In the Ji = 0 limit, the model becomes a 



FIG. 4: Illustration of the model Eq. (141 in terms of Majo¬ 
rana fermions. After the Jordan-Wigner transform, the com¬ 
plex fermion Ci on each site (yellow block) can be split into 
two Majorana fermions (white circles) as Ci = Xii+iXi 2 . Then 
the model becomes two decoupled free Majorana chains. 












7 


classical Ising model with random Zeeman field, 

H = -Y^ (15) 

i 


which is already solved in the Ising basis. However this 
model can not be mapped to a free fermion model, so the 
interaction effect still remains (although rather trivial). 
The model Eq. (151 has only one phase (the PM phase), 
where the strong-if and the strong-/i limits are smoothly 
connected without phase transition. 


B. Characterization of the RG Flow 


As an example, we will take the quantum Ising model 
in Eq. (11) (or equivalently the interacting fermion model 
in Eq. (12)), and consider the following initial distribu¬ 
tions of the coefficients Ji, Ki and hi, 


J \ 


l/Fo 


P{J)dJ = — \^—j dJ for J e [0, Jo], 

1 / K \ 

P(if)dX= — f—j dKioTKG[0,Ko], (16) 


P{h)dh = 


1 

Pq/i V h, 


'-0 


l/Fo 


dh for h € [0, ho]. 


The disorder strength is controlled by a single parameter 
Pq. Larger Pq corresponds to stronger disorder. These 
initial distributions of the coefficients are expected to flow 
under SBRG. 

The first question to ask is whether the distributions 
flows towards the strong disorder limit? The answer how¬ 
ever depends on the phase of the Hamiltonian. If the 
Hamiltonian is in the MBL phase, SBRG will flow to¬ 
wards strong disorder. On the other hand, if the Hamil¬ 
tonian is in the ETH phase, SBRG will flow away from 
the strong disorder limit. 

To quantify the RG flow, we investigate the many- 
body Thouless parameter introduced in Ref.[89H9TJ The 
Thouless parameter can be viewed as a many-body gen¬ 
eralization of the Thouless conductance^ in the single¬ 
particle Anderson problem, which is the ratio of the off- 
diagonal resonance energy Vmn to the diagonal level spac¬ 
ing, as g = \V,nn\/\Em - En\, where Kn„ denotes the 
matrix element of a local perturbation represented in 
the many-body eigen basis. To give a crude estimate 
of the Thouless parameter in the SBRG implementa¬ 
tion, we take out the Hamiltonian Eq. @ at each RG 
step: H = Hq -|- A -|- E, which contains the leading en¬ 
ergy scale term Hq and many block-off-diagonal terms 
E = El -I- E 2 . For each block-off-diagonal term, we 
collect the ratio of its energy scale to the leading energy 
scale, and define the ratio as the Thouless parameter 


a = 


lis.ll 

\\HoW 


(17) 


This ratio lies between zero and one, and is the small pa¬ 
rameter to control the perturbative treatment in the RG 
scheme. In the strong disorder limit, the leading energy 
scale is expected to be much larger than all the other 
terms in the Hamiltonian, and g tends to zero. Physi¬ 
cally this ratio also resembles the ratio of the resonance 
energy scale over the many-body level spacing. At each 
RG step, the spectrum bifurcates. The leading energy 
scale ||i7o|| is the amount of spectrum splitting at the 
current RG step, which also controls the level spacing. 
The block-off-diagonal terms E^ flips the emergent qubit 
at the current RG step, which can be used to characterize 
the resonance energy scale. 



/ = log2 {NoIN) 

FIG. 5: Distribution of the logarithmic Thouless parame¬ 
ter G (plotted in logarithmic scale) for disordered quantum 
Ising model at (a) {Jo,Ko,ho) = (1,1,1) and Fq = 1, (b) 
{Jq, Ko,ho) = (2,1,1) and Fq = 1. The calculation is per¬ 
formed on a 512-site lattice over sufficient amount of random 
realizations. G is collected when the RG runs to the steps 
that the number of physical qubits is reduced to N. (c) Flow 
of the typical G under RG. {Nq/N) characterizes the length 
scale. Smaller N corresponds to longer length scale. 


It will be more convenient to study the logarith¬ 
mic Thouless parameter, denoted as G = Ing = 
ln(||Ei||/||iJo||). We collect the parameter G with the 
RG flow for various disorder realizations, from which we 
can study how the probability distribution P{G)dG flows 
under RG. Fig.j^a) shows the case that the system is in 
the MBL phase. The distribution P{G) ~ has an 
exponential tail in the large negative G limit (G —)■ — 00 ). 
The tail keeps broadening with the RG flow (correspond¬ 
ing to P flowing towards -l-oo). Therefore the average G 
will flow to — 00 , or the typical value of the Thouless pa¬ 
rameter g will flow to zero, as shown in Fig. I^c). This 
demonstrates that SBRG flows towards the strong dis¬ 
order limit in the MBL phase, and become progressively 














accurate. In contrast, Fig.j^b) shows the case that the 
system is close to the ETH phase (or weakly thermal- 
ized). In this case, the distribution of G almost does not 
flow under RG, and even shifts (to the right) towards 
G —?► 0 near the end of the RG flow. Gorrespondingly av¬ 
erage G does not flow towards (or even flows away from) 
the strong disorder limit, as shown in Fig. I^c). So the 
RG will eventually break down in the ETH phase, and 
the signature of thermalization can be indicated from the 
trend of the RG flow. 

In conclusion, SBRG only works for MBL systems, 
and can not approach the thermalization transition or 
enter the ETH phase. Whether SBRG applies to the 
marginal MBL system is more subtle. Marginal MBL 
states are more delocalized than MBL states, and hence 
easier to thermalize if the disorder strength is weaken. It 
is still an open question whether marginal MBL states 
are stable against thermalization in the presence of in¬ 
teraction. However even if the marginal MBL state is 
unstable against thermalization, the thermalization ef¬ 
fect is expected to be weak in the strong disorder regime. 
So it is still valid to talk about marginal MBL states as 
the short-range physics (or for finite-sized system). In 
this sense, SBRG is applicable to the marginal MBL sys¬ 
tem as well. Various universal scaling properties of the 
marginal MBL state can be studied by SBRG, as will be 
shown in later sections. 


C. Benchmarking the Many-Body Spectrum 


At the end of the RG flow, we will obtain the full en¬ 
ergy spectrum, encoded in the effective Hamiltonian 

HeS = '^eiTi+'^eijTiTj+ ^ eijkTiTjTk-\ -, (18) 

i i<j i<j<k 


where = ±1 is the emergent conserved quantities iden¬ 
tified in SBRG flow. By iterating over all configurations 
of Ti, all the eigen energies can be obtained from Eq. (18) 
approximately. 



Eed 


FIG. 6: Benchmark the many-body spectrum obtained by 
SBRG with ED. Galculation performed on a 8-site lattice, 
with Ji, Ki, hi independently drawn from the distribution in 
Eq. (161 with Jo = Kq = ho — 1 for Fq = 1, 2, 5. 


The many-body energy spectrum estimated by SBRG 
can be benchmarked with ED on small-sized system (on 
which ED can be performed). In Fig.[^ we show such 
a comparison of the many-body spectrum in an 8-qubit 
system (256 eigen energies) for various shapes of the ini¬ 
tial distributions of (controlled by Fq = 1,2, 5, where 
the larger Fq the stronger disorder). In general, all points 
almost line up straightly along the diagonal line, showing 
good agreement between SBRG estimation and the ex¬ 
act spectrum from ED. Also as the disorder gets stronger 
(yellow points Fq = 5), the agreement becomes better, 
showing that the accuracy of SBRG is systematically im¬ 
proved towards the strong disorder limit. 


D. Statistics of Energy Coefficients 


Having obtained the effective Hamiltonian in Eq. (18) 


from SBRG, we can look at the statistics of the energy 
coefficients e^, e^, e^-fc ect. Let us denote the n-body 
coefficient ^(n)> study its probability dis¬ 

tribution P(e(„))de(„). The numerical results are shown 
in Fig.0 


(Jo,Ko, bo):(1.1.0) (1,1,1) (2,1,1) 



FIG. 7: Probability distribution of the n-body energy coef¬ 
ficients £(„) in the effective Hamiltonian. The statistics is 
collected from 200 random realizations on a 128-site lattice 
for each case. Different columns correspond to different set¬ 
tings of the initial scales {Jo, Kq, ho), and different rows cor¬ 
respond to different randomness of the initial distributions. 
For {Jo, Ko,ho) = (1,1,0), the many-body interaction terms 
C{n) (n = 2,3,4, • ■ •) are not generated, and hence their statis¬ 
tics are not shown. Fq = 1 corresponds to the uniform initial 
distribution, and larger Fq corresponds to stronger disorder. 


First let us look at the (Jo,iFoi^o) = (IjliO) free 
fermion critical point, mentioned in Eq. (|14[), as shown 
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in the left column of Fig.[^ In this case, the many- 
body terms e(„) (n = 2,3,4, •••) are not generated un¬ 
der the RG flow, showing that SBRG has automatically 
found the single particle basis to diagonalize the Hamil¬ 
tonian. The distributions of the single particle levels e(i) 
exhibit the Dyson singularitj^sHIlll Jqw energy (or low 
frequency), i.e. P(e(i)) ~ eC^^(ln(IT/e(i)))“^ where W is 
the typical band width. 

The middle column of Fig.[^ shows the cases of 
(Jo, Kg, ho) = (1,1,1), which is in the PM (or SPT triv¬ 
ial) MBL phase according to phase diagram Fig.[^ (to 
be discussed later). For the weak disorder case Fq = 1, 
the single-particle DOS (the distribution of e(i)) drops 
to zero at the low frequency limit (e(i) —>■ 0), which is in 
reminiscence of the fermion single-particle band gap in 
the clean limit. As the disorder gets stronger (Fq = 2, 5), 
the band gap is gradually filled up by the in-gap local¬ 
ized states. The right column of Fig.[^ shows the cases 
of {Jo,Ko,ho) = (2,1,1), which sits right at the inter¬ 
acting marginal MBL critical line. The single-particle 
DOS remains gapless in these cases, which is consistent 
with the fact that the single-particle gap must close at 
the transition between the topological (zz = 1) and the 
trivial (v = 0) ID fermion SPT phaseseven in the 
presence of disorder and interaction. 

Now let us turn to the many-body terms £(„) (n = 
2, 3,4 • • •). Away from the free fermion limit, the many- 
body terms will appear under the RG flow, and the ef¬ 
fective Hamiltonian iJeff becomes an interacting one, as 
shown in the middle and right columns of Fig.[^ In 
the strong disorder limit Pq = 5, the distributions look 
similar to the free fermion case, and little interaction 
is generated among the emergent conserved quantities. 
As the disorder gets weaker, the interaction becomes 
stronger (as the distribution shifts towards higher energy 
scale) and involves more emergent conserved quantities. 
Many of the interaction terms are UV-IR mixing terms, 
which describe how the spectrum branching in the (high- 
frequency) UV limit could shift and rearrange the energy 
levels in (low-frequency) IR limit. 

The many-body energy coefficients €(„) also set the en¬ 
ergy scale of many-body resonances in the MBL system. 
To quantify the n-body resonance energy scale, we stud¬ 
ied the norm of the n-body terms ||e(n)||, defined as 



We found that the n-body resonance energy scale decays 
with n exponentially, 

||e(n)||-^t^e-/«, (20) 

as shown in Fig.lM This observation justifies the same 
proposal in Ref.|^ From the exponential decay of the n- 
body resonance energy scale, it was further arguecpSl that 
the conductivity should follow a power-law of frecmency 
(t(w) ~ in the deep MBL phase. From Fig.[^ one 


can also observe that as the disorder gets weaker (smaller 
Fq), ||e(„) II decays slower with n, meaning that the many- 
body interaction/resonance will become more important 
for weaker disorder. 


(Jo,/<o, ho) = (1,1,1) 



n 

FIG. 8: The energy scale ||e(„)|| of n-body terms decays with 
n exponentially in the MBL phase. The decay rate depends 
on the initial randomness Fq. 


In the strong disorder and weak interaction limit, 
where the many-body terms are suppressed, the effective 
Hamiltonian JLeff = i® expected to produce the 

Poisson level statistics due to the lack of level repulsion, 
which is consistent with the MBL physics. Away from 
that limit, the many-body terms can reshuffle the en¬ 
ergy levels at low-frequency, which may change the level 
statistics drastically. However, it is not clear to us if the 
level statistics imposes any constraint on the distribu¬ 
tions P(e(„)) of the energy coefficients, or if it is possible 
to determine the level statistics from P(e(„)). The dis¬ 
cussion along this line may go beyond the scope of this 
work, and we will leave these interesting questions for 
future study. 

E. RG Transform as Clifford Circuit 

SBRG is a Hilbert space preserving RG. The RG trans¬ 
formation contains no isometry, but only unitary trans¬ 
forms. So by collecting the unitary transforms that have 
been performed in each RG step, we can obtain a uni¬ 
tary mapping which maps the physical Hilbert space to 
the emergent Hilbert space, and consequently maps the 
initial Hamiltonian H to the effective Hamiltonian Hgg 
(approximately). In each RG step, the Hamiltonian is 
conjugated by two unitary transforms: a Glifford group 
rotation Rk followed by a Schrieffer-Wolff transformation 
Sk (here fc = I, 2, • • • labels the fcth RG step). Therefore 
the whole RG transformation is a product of them as 

I/rg = R1S1R2S2 • • • = n ( 21 ) 

k 

which approximately diagonalize the many-body Hamil¬ 
tonian H Hes — U^QHU^^G. Note that in SBRG 
algorithm, the Schrieffer-Wolff transformation Sk is not 
performed exactly, but only carried out to the 2nd order 
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in the perturbative expansion, where all the higher order 
terms are truncated. So the truncation error can build 
up in the Hamiltonian under the RG flow, and hence 
the resulting RG transform [/rg is only an approximate 
unitary transform to diagonalize the Hamiltonian. Nev¬ 
ertheless the approximation is expected to be good if the 
disorder is strong (such as in the full MBL phase). 

Therefore the RG transform C/rg actually encodes all 
the many-body eigenstates of H to some approximation. 
To retrieve the eigenstates from C/rg, we just need to take 
the eigenstates of the effective Hamiltonian Hgg, which 
are direct product states |{ri}) = |ti) (g) |t 2 ) g ••• of 
the emergent qubits ti = ±1, and transform them back 
to the original basis by reversing the RG transformation 
~ C/RG|{ri}). In this sense, C/rg can be viewed as 
a quantum circuit!^ which prepares the entangled eigen¬ 
state from the disentangled product state. Different en¬ 
ergy eigenstates in the MBL spectrum correspond 

to different input product state |{ri}) of the emergent 
qubits. 

Although we have the Clifford rotation Rk and SW 
transformation Sk matrices in hand, computing their 
product explicitly as in Eq. (21) is still expensive in nu¬ 


merics. So at this stage, we make further approximations 
by dropping all the SW transformations Sk, and recon¬ 
struct the eigenstates just from the Clifford rotations Rk, 
because Sk are expected to be close to identity. With Rk 
only, the remaining unitary transformation is a quantum 
circuit of random Clifford gates, dubbed as th e rando m 
Clifford circuit, or the random stabilizer circuitP^^^^^^ 


Uci = Rk- 


( 22 ) 


Each Clifford gate Rk here is a controlled gate in gen¬ 
eral, which transforms one physical qubit to one emergent 
qubit under the control of the previously identified emer¬ 
gent qubits, as illustrated in Fig.[^ The emergent qubits 
will not be further rotated by the later gates, but may 
serve as the control qubits for the later gates. The Clif¬ 
ford circuit is a further approximation of the RG trans¬ 
formation, which gives cruder estimate of the eigenstate 

\'i'ir.})^Uci\{T.}). (23) 

The accuracy is traded for efficiency. Since Uci\{Ti}) is 
a stabilizer state, many physical properties can be effi¬ 
ciently calculated using the stabilizer formalism. 

It worth mention that dropping the Schrieffer-Wolff 
transformations Sk only affects the accuracy of the esti¬ 
mated eigenstates, but not the SBRG flow, because the 
RG flow is purely guided by the flow of the Hamilto¬ 
nian without referring to Sk- In particular, the 2nd or¬ 
der perturbation can be performed on the Hamiltonian 
level according to Eq. without explicitly calculating 
the Schrieffer-Wolff transformation S in Eq. ([^. So we 
can first determine the Hamiltonian flow under SBRG, 
and then reconstruct the Clifford circuit C7ci indepen¬ 
dently after the RG flow. 


emergent quibits 



FIG. 9: The random Clifford circuit Uci- Each Clifford ro¬ 
tation Rk is viewed as a (controlled) quantum gate, where 
the black dots mark the control qubits and the yellow squares 
are generic unitary gates. Under the action of each gate, one 
physical qubit (in blue) will be transformed to one emergent 
qubit (in red). Only the emergent qubit serves as the control. 


The Clifford circuit actually provides the MRS rep¬ 
resentation for all eigenstates. Because each Clifford 
gate can be represented as a matrix product operator 
(MPOJP^with fixed bond dimension 2. Then the entire 
Clifford circuit is a large MPO with total bound dimen¬ 
sion bounded by the circuit depth. As a RG transforma¬ 
tion, the circuit depth is typically at most logarithmic 
in system size. The logarithmic depth is only saturated 
for the marginal MBL system. For systems deep in the 
MBL phase, the RG will terminate at finite depth, and 
the resulting MPO actually has a (smaller) bound di¬ 
mension independent of the system size.^ Inputting a 
direct product state from the emergent qubit side corre¬ 
sponds to fixing the emergent legs of the MPO (red legs in 
Fig.§, which turns the MPO to an MPS representation 
of the MBL eigenstate. Thus once the random Clifford 
circuit is generated from the SBRG flow, we have also ob¬ 
tained the (approximate) MPS representation for every 
energy eigenstate in the whole spectrum. Because the 
bipartite entanglement entropy is bounded by the bound 
dimension, the entanglement entropy can scale logarith¬ 
mically at margi nal MBL, and follows area-law in the 
MBL phase.li^ 


F. Benchmarking the Many-Body States 


To check how good is the approximation of the ran¬ 
dom Clifford circuit, we can benchmark the circuit with 
the exact diagonalization (ED) unitary transformation. 
We take the quantum Ising model in Eq. 0 on a 8- 
site lattice. First we run SBRG and collect the Clifford 
rotations Rk {k = 1, - - - ,8). Then we assemble the Clif¬ 
ford circuit Uci = nfe=i Rk and retrieve the eigenstates 
from Eq. (231. As we enumerate over all configurations of 


the emergent conserved quantities {ri}, we can obtain all 
the approximate eigenstates Finally we overlap 

these states with the exact energy eigenstates obtained 
by ED. Typically one can order the eigenstates by their 
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eigen energies. But since the energy spectrum estimated 
by SBRG also contains small error, a few eigenstates 
might be reshuffled. Therefore we match the correspond¬ 
ing eigenstates by the maximal overlap (i mple mented by 
the maximal bipartite matching algorithrrfl®) , and cal¬ 
culate the fidelity f^n} = where 

is the ED eigenstate that matches. The fidelity of each 
RG-generated eigenstate is plotted against its energy in 
Fig.[^ For moderate randomness Tq = 1, the fidelity 
is typically around 0.5, and improves systematically as 
the randomness gets stronger (see the Tg = 5 data). It 
is fair to say the result is acceptable for the wave func¬ 
tion overlap in the 256-dimensional many-body Hilbert 
space. The fidelity level of SBRG is also comparable 
to that of RSRG-X (as benchmarked in Ref. ElTIl . show¬ 
ing that SBRG has a similar performance as RSRG-X. 
Admittedly, SBRG is a fast but less accurate method. 
Other approaches such as DMRG-XP^^^ or tensor net¬ 
work optimizatiorPSm can achieve much higher eigen¬ 
state fidelity. If high-quality eigenstates are desired, the 
eigenstate estimated by SBRG can be passed as the initial 
MRS state to DMRG-X for further refinement. We also 
observe from Fig.that the fidelity is typically higher 
for the ground state and top state in the spectrum, and 
lower for the states in the middle of the spectrum. This is 
because the middle state has a stronger tendency to ther- 
malize, and we know that SBRG and the Clifford circuit 
will become inaccurate towards thermalization. In con¬ 
clusion, as long as we are in the MBL (or marginal MBL) 
phase, the Clifford circuit can provide reasonable approx¬ 
imation for the many-body eigenstates in the whole spec¬ 
trum. 



^RG 


FIG. 10: Many-body fidelity of the energy eigenstates ob¬ 
tained by the Clifford circuit benchmarked with ED. Cal¬ 
culation performed on a 8-site lattice, with Ji, Ki, hi in¬ 
dependently drawn from the distribution in Eq. (161 with 
Jo = Ao = ho = 1 for Eo = 1, 2, 5. 


For larger-sized system where the full-spectrum ED is 
no longer available, we can check the accuracy of the 
Clifford circuit by the energy variance, as proposed in 
Ref.[5H We define the energy variance averaged over the 


approximate eigenstates as 

(24) 

where H is the initial model Hamiltonian. If were 

exact eigenstates, the mean energy variance SE"^ should 
vanish. So the non-vanishing SE^ indicates how far the 
Clifford circuit deviates from the exact diagonalization. 
However, SE"^ grows with the system size because it has 
the dimension of the total energy squared and the to¬ 
tal energy is an extensive quantity that scales with the 
system size. To get rid of t he system size dependence, 
we choose to normalize SE"^ by the total variance of the 
energy 

( 25 ) 

and the result is shown in Fig. El For various model pa¬ 
rameters (JojA'ojho), the normalized mean energy vari¬ 
ance 6E'^/E^ is controlled below ~ 0.4,hJ^ and decreases 
all the way to zero towards the strong disorder limit 
where the Clifford circuit becomes exact. 



FIG. 11: Normalized energy variance SE'^jE^ v.s. the ran¬ 
domness Fo of the initial distribution, for various parameters 
(Jo, Ao, ho). The calculation is performed on a 64-site lattice 
with 100 random realizations for each data point. 


In conclusion, SBRG generates two important sets of 
data: a fixed-point Hamiltonian i/gS which gives the full 
energy spectrum to the 2nd order perturbation, and a 
Clifford circuit Uci which encodes all the eigenstates to 
the 0th order perturbation. Admittedly, the Clifford cir¬ 
cuit may be a crude approximation in terms of encoding 
the many-body eigenstates, but it has nice properties to 
enable highly efficient calculation of many physical quan¬ 
tities, which will be discussed in the following section. 


IV. ENTANGLEMENT HOLOGRAPHIC 
MAPPING 

A. SBRG and Holography 

SBRG also provides a holographic interpretation of 
the MBL system. This follows from the idea that every 
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Hilbert-space-preserving RG can be interpreted as a holo¬ 
graphic mapping, and the RG transformation is a unitary 
mapping from the holographic boundary to holographic 
bulk. The bulk qubits defined by this mapping can be 
considered as local degrees of freedom in an emergent 
bulk geometry, whose radial direction (or the holographic 
extra dimension) corresponds to the energy (frequency) 
scale, running from UV to IR, as illustrated in Fig.[T^ 
The close relation between RG and holography has been 
actively studied in both condensed matter ph ysics and 
high energy physics communities .ISMl o i io il tlMI SBRG is 
also a Hilbert-space-preserving RG, so it must also have 
a holographic interpretation. 

Here we will follow the construction of EHM proposed 
by one of the authors in Ref. 18118^ The construction 
was first applied to a free fermion system without disor¬ 
der. The model describes spinless fermions hopping on a 
ID lattice. At each RG step, two neighboring sites are 
coarse grained to an effective site either in the low-energy 
subspace or in the high-energy subspace, depending on 
the state of the bulk fermion (either empty or occupied). 
The coarse graining is implemented as unitary transfor¬ 
mations of the fermion basis, represented as small yellow 
disks in Fig.[T^a). In this way, a tensor network is gen¬ 
erated under RG, which maps the boundary fermions to 
the bulk. 

Now for SBRG, RG steps are implemented by Clif¬ 
ford gates (approximately), depicted as yellow blocks 
in Fig.[^b). Each Clifford gate identifies an emergent 
qubit on one of its output leg. We pull out the emer¬ 
gent qubits as bulk degrees of freedom, because they are 
the ones who control the spectrum branching at each 
RG step. The bulk qubits reside at the radius position 
determined by their energy scales at which they are iden¬ 
tified by SBRG. The energy scale can be easily read off 
from the effective Hamiltonian Heff = + ’ ’ ’ ? such 

that the energy scale associated to n is simply a, i.e. 
the single-body energy coefficients in iJeff- In this way, 
the random Clifford circuit is interpreted as the EHM 
network. Since we are dealing with disordered systems, 
our EHM network is not a regular tensor network, and 
the bulk geometry also inherits the randomness on the 
boundary. The tensor network maps the MBL eigen¬ 
states on the boundary to a direct product state in the 
bulk, so it is also a disentangler network or a random 
version of MERA,1^ which progressively removes the en¬ 
tanglements in the wave function. 

An important feature of the SBRG-generated EHM 
network is that the transformations performed in the IR 
region are controlled by the states of the emergent qubits 
in the UV region, as depicted in Fig.[I^b) by the gray 
lines connecting the IR gates to the UV emergent qubits. 
The emergent qubit state in the UV region determines 
the large-scale spectrum branching and roughly locates 
the system around certain energy density in the spec¬ 
trum, therefore the IR transformations can be affected 
by the choice of the qubit state that has been made in the 
UV layer. Correspondingly, the spectrum branching de- 



FIG. 12: EHM network of (a) free fermion system without 
disorder,^ and (b) MBL system. Both EHM networks map 
the physical qubits (in blue) on the holographic boundary to 
the emergent qubits (in red) in the holographic bulk. The 
unitary transform in each RG step is represented by a yellow 
block. Each RG step will identify a new bulk qubit. 


pendence of the IR physics is reflected in the many-body 

terms -b H-in the effective 

Hamiltonian as a UV-IR mixing effect. 

The idea of holography is useful, because it provides 
a geometric interpretation of the entanglement struc¬ 
tures in the quantum many-body state. For exam¬ 
ple, the entanglement entropy can be interpreted as 
the minimal surface in the holo graphic bulk following 
the Ryu-Takayanagi formula.l^^ The correlation func¬ 
tion or mutual information is also related to the holo¬ 
graphic geodesic distance. In general, a full-spectrum 
holographic mapping for generic many-body system is 
challenging. However the MBL systems are special, in 
the sense that they are “quasi-solvable”, which allows 
a Hilbert-space-preserving RG and a controlled holo¬ 
graphic mapping of the entire many-body Hilbert space. 
In the following we will discuss several physical proper¬ 
ties of MBL systems, and the holographic geometry in 
the bulk. 


B. Stablizer Properties 


Let us first look at the emergent conserved quantities 
Ti, also known as th e localized bit^^ or the local inte¬ 
grals of motiorpilSSl^ which play important roles in the 
phenomenology of MBL systems. On one hand, they la¬ 
bel the emergent qubit states in the holographic bulk. On 
the other hand, they are the stabilizers (commuting pro¬ 
jectors) for the eigenstates on the holographic boundary. 
The representation of the emergent conserved quantities 
in the physical Hilbert space fi can be approximately 
found by inversely applying the Glifford circuit, 

= (26) 


where [3^] denotes the set of Pauli indices of the form 
[0 • • • ]3[0 • • • ] with 3 appears at the Ith qubit position for 
z = 1,2, • • • ,N. We use a hat to emphasize that the oper¬ 
ator Ti is represented in the physical Hilbert space, acting 
on the holographic boundary. From the effective Hamil¬ 
tonian Hgff in Eq. (181, we know that every eigenstate in 
the physical Hilbert space is a stabilizer state (approx- 
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imately), stabilized by Ti to the value of the emergent 
conserved quantity = ±1 


Vi : =Ti\'i’[r,})- 


(27) 


Each stabilizer fi is also a Pauli operator, because the 
Clifford circuit can only take the Pauli operator to 
another Pauli operator. To gain some intuition of the 
stabilizers, let us plot them in Fig.[^ We start with 
the quantum Ising model in Eq. (11) (or as the fermion 
model in Eq. (|l2|)), and take the random distribution 
from Eq. (16) with the parameters ( Jq, Kq, ho) = (1,0,1) 
for Fig.|13[a) and (Jo,^o,^o) = (1,1,1) for Fig.j^b). 
The stabilizers are collected along SBRG flow, and re¬ 
stored to the physical Hilbert space by Eq. (26 1 . 
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FIG. 13: Snapshots of the stabilizers ti arranged vertically by 
their energy scales ei from UV (bottom) to IR (top), for (a) 
{Jo,Ko,ho) = (1,0,1) and (b) {Jo,Ko,ho) = (1,1,1). Hor¬ 
izontal axis is the real-space position in a 32-site subsystem 
taken from the 256-site lattice. Each stabilizer is represented 
as a product of Pauli matrices (grouped in gray shadows), 
and each Pauli matrix is illustrated by a color dot (red= a^, 
green= a^, yellow^ a^). 


One important question is the locality of these sta¬ 
bilizers. As we can see from Fig. |13[ most of them are 
quite well localized. Fig.[T^a) shows the scenario at the 
marginal MBL criticality (at the critical point of free Ma- 
jorana/Ising chain). In the UV regime (i.e. — In \ei\ ^ 0), 
many stabilizers are of the form (two red 

dots) or tT [0 "] 3 [o-"] single yellow dot), corresponding 
to the strong bond or the strong site operators respec¬ 
tively. As the RG flows towards the IR limit (larger 
— ln|ei|, lower frequency), longer stabilizers are found, 
corresponding to the longer-range coupling of the Majo- 
rana fermions. The fermion Jordan-Wigner string can 
also be observed as crI^ -1 (a line of yellow dots between 
the red/green dots). These long-range couplings are quite 
rare and have very low energy scale, which follow the uni¬ 
versal scaling behavior — Inje^l ^ with li being the 
length of the stabilizer.l^ We can tune the system away 


from the criticality to the MBL phase by the fermion 
interaction Ki{ni — |)(ni_|_i — 5 ). Because on the mean- 
field level, the interaction effectively enhances the ran¬ 
domness of the site field as hi ^ hi + Ki{ni+i — |), 
and therefore tips the balance between the bond Ji 
and the site hi terms.^ As shown in Fig. 13 b), with 


{ Jo, Ko,ho) = (1,1,1), all the stabilizers become short- 
range and concentrated near the UV regime. We see that 
SBRG stops flowing within a few orders of the energy 
scale in the MBL phase. 


(Jo, Ko, ho) 
— A(1,0,1) 

— B(1,1,0) 
— D(4,1,1) 

— E(2,1,1) 

— F(1,1,1) 


0 10 20 30 40 

/ 




FIG. 14: Probability distribution P{1) of the stabilizer length 
I, shown as (a) log-linear plot and as (b) log-log plot. The 
statistics is collected from 1500 random realizations on a 128- 
site lattice for each par am eter {Jo, Ko, ho) as marked out in 
the phase diagram Fig. 17 accordingly, with Pq = 1. In (b), 
the lines are shifted vertically from each other (by 4 vertical 
units) for clarity. 


Having the above picture in mind, we now check the 
stabilizer locality quantitatively. First, we collect the 
stabilizers obtained by SBRG, and count the number of 
stabilizers of each length 1. The length of the stabilizer 
is defined as the real-space separation from the left-most 
qubit to the right-most qubit that are acted by the sta¬ 
bilizer in the physical Hilbert space. For each initial 
scales {Jo,Ko,ho) of the couplings, we repeat the cal¬ 
culation for many random realizations to obtain reliable 
statistics. The statistical frequency is then normalized 
to the probability distribution P{l)dl of the stabilizer 
length 1. Fig. 1^ a) shows the probability distribution 
P{1) v.s. I in the log-linear plot. According to the phase 
diagram Fig. 17 to be shown latter, D{A, 1,1) belongs to 
the SG MBL phase, and F(l,l,l) belongs to the PM 
MBL phase. We can see, in both MBL phases, the prob¬ 
ability distribution decays exponentially P{1) ~ 
meaning that the stabilizers are exponentially localized 
in the MBL phases. Our result can be viewed as an sup- 
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portive evidence for Ref. l21l2^ and also agrees with the 
analysis in Ref.[28l 

However at the marginal MBL critical point, the sta¬ 
bilizers are no longer short-ranged, but become quasi- 
long-ranged with the power law distribution P{1) ~ l~°‘. 
According to the phase diagram in Fig.[^ A(1,0,1), 
H( 1 , 1 , 0 ) and A( 2 , 1 , 1 ) are on (or very close to) the 
marginal MBL critical line. Then in Fig. [Id]) a), their 
probability distributions P{1) indeed deviate from the 
exponential decay, and develop long tails for large 1. 
To confirm the power law behavior, we switch to the 
log-log plot in Fig.[l4](b), and find that the data fol¬ 
low straight lines nicely. For the free fermion chains 
A( 1 , 0 , 1 ) and i3(l,l,0), the universal exponent a = 2 
is expected, because according to RSRG fixed point 
solution,!^ the length scale I is related to the logarith¬ 
mic energy scale F = — Ine as / ^ F^, so the number 
of stabilizers identified at that energy scale should fol¬ 
low P{l)dl = {No/l){dr/T) ^ l~^dl where Nq is the ini¬ 
tial number of sites (qubits). We confirm the exponent 
a = 2 in Fig. mb) (see the lines A and B). For the 
interacting fermion chain E{2, 1,1) at criticality, the RG 
fix-point theory is not known to us. Our numerical result 
seems to imply roughly the same exponent a ~ 2.1 as the 
line E in Fig. mb)- The numerical calculation obtains a 
slightly larger exponent, because it is difficult to exactly 
locate the critical line, thus the point E{2, 1,1) may be 
a little off-critical, therefore P{1) is subject to a weak 
exponential decay, which could result in a seemly larger 
exponent from the power-law fitting. In Tab.|l| we col¬ 
lect the exponents a calculated on a few points along (or 
near) the marginal MBL critical line. We do not observe 
systematic or significant deviation from a = 2. So our 
result implies that the lengths of the stabilizers are likely 
to follow the same universal scaling behavior P{1) ~ 
along the whole critical line. 


TABLE I: The scaling exponent a obtained by the power-law 
fitting along (or near) the phase boundary. 


{Jo,Ko, ho) 

a 

(Jo, Kq, ho] 

a 

(1,0,1) 

2.03 ±0.09 

(5,3,2) 

2.13 ±0.03 

(8,1,7) 

2.03 ±0.09 

(6,4,1) 

2.12 ±0.06 

(4,1,3) 

1.98 ±0.09 

(7,5,1) 

2.06 ±0.06 

(2,1,1) 

2.11 ±0.07 

(1,1,0) 

1.99 ±0.06 


Another important question is the dynamical scaling at 
the marginal MBL criticality. The absence of diffusion 
in the MBL system is characterized by the infinite dy¬ 
namical exponent z —> oo in ^ ~ . Or more precisely, 

this implies the logarithmic dynamical scaling I ~ (Int)** 
at the marginal MBL critical point, where 77 = 2 is ex¬ 
pected from the previous RSRG study on the transverse 
field Ising model.l^ We can check the dynamical scal¬ 
ing by studying the energy scale v.s. the length scale of 
the stabilizers identified by SBRG. For every fi, we 
read off its energy scale Ci from the effective Hamilto- 



(Jo, Ko, ho) 

• A(1,0,1) 

• B(1,1,0) 

• E(2,1,1) 


FIG. 15: Universal dynamical scaling. Different colors repre¬ 
sent the different settings of the initial scales {Jo, Kq, ho) as 
marked out in the phase diagram Fig.[^accordingly. The cal¬ 
culation is performed on a 128-site lattice with 1000 random 
realizations. 


nian e^Ti and measure its length U in the 

physical space. We collect the pair (eiji) for all stabi¬ 
lizers in the system, and repeat the calculation for many 
random realizations to obtain sufficient amount of sam¬ 
ples. Fig.[^ shows that the data points roughly follow 
the — Ine ~ Vl behavior, and hence the dynamical scal¬ 
ing I ^ (—Ine)^ ^ (Int)^ is justified. Our data imply 
that the universal dynamical scaling not only applies to 
the free fermion chains A(l, 0,1), B{1, 1,0), but also ap¬ 
plies to interacting fermion chain E{2,1,1). Thus we 
conclude that for the quantum Ising model, the whole 
marginal MBL critical line is characterized by the same 
universal dynamical scaling with z —^ 00 . 


C. Entanglement Entropy 


In the Glifford circuit formalism, every energy eigen¬ 
state in the physical Hilbert space is approximated by a 
stabilizer state as in Eq. (27) that corresponds to 

the direct product state 11^}) in the holographic bulk. In 
this case, all the entanglement structure of the stabilizer 
state is given by the Glifford circuit, and the bulk qubits 
will make no contribution because they are completely 
disentangled. However, we know the Glifford circuit only 
provides the 0th order approximation of the eigenstates. 
If we fix the Clifford circuit, and map the exact physical 
eigenstate from the holographic boundary into the bulk, 
then the corresponding bulk state will actually has some 
entanglement. The entangled bulk qubits will provide an 
additional layer of entanglement on top of the entangle¬ 
ment given by the Clifford circuit. Therefore the EE Se 
of an energy eigenstate is given by the background en¬ 
tropy ^ci based on the Clifford circuit plus the correction 
S'buik due to the bulk qubit entanglement (which can be 
negative), 


Se = S'ci + -S'buik- (28) 

A nice property of the Clifford gates Rk is that they are 
all of the same MPO bond dimension Du = 2. So the 
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MPO bond dimension of the Clifford circuit is simply 
D = = 2‘^ with d being the circuit depth. This im¬ 

poses an upper bound on the EE of the Clifford circuit: 
S'ci < InZl = dln2. In the MBL phase, the depth of 
the Clifford circuit is finite, and S'ci follows the area-law 
scaling. At the marginal MBL criticality, the depth of 
the Clifford circuit grows with the system size logarith¬ 
mically, so Sci ~ InL also follows the logarithmic scal¬ 
ing. Typically the Clifford circuit generated by SBRG 
can not provide a volume law entanglement scaling, be¬ 
cause for an A^-qubit system, there are only N stabilizers, 
and most of them are localized. Therefore the volume law 
entanglement (if any) can only come from the holographic 
bulk as Sbuik ~ L. But in that case, Sci is overwhelmed 
by Sbuik SO that the Clifford circuit is no longer a good 
starting point to encode the volume-law state. 

In the strong disorder limit, the bulk contribution Sbuik 
can be neglected, and the EE can be estimated from the 
Clifford circuit alone, i.e. Se — Sci- It turns out that 
the biparti te E E can be calculated efficiently for stabi¬ 
lizer states,^^^ and the method is reviewed in Appendix 
[B| For each Clifford circuit generated by SBRG, we can 
choose a subsystem of the length L and calculate the EE 
associated to its reduced density matrix using the stabi¬ 
lizer formalism. Then we repeat the calculation with the 
entanglement cuts translated through out the system to 
average over the disorder configuration. We also repeat 
for many random realizations of the system to get more 
reliable statistics. 


(Jo, Ko, ho) 

• A(1,0,1) 

• B(1,1,0) 

• C(1,0,0) 
. E(2,1,1) 
. F(1,1,1) 

• G(0,1,1) 

0 1 2 3 4 5 6 7 

lOQa L 



FIG. 16: Scaling of EE Se for various initial (Jo, Kq, ho) with 
Fo = 1. L is the length of the subsystem for the EE measure¬ 
ment. The points are marked out in the phase diagram Eig.|17| 
correspondingly. The shading denotes the confidence interval 
of one standard deviation. The calculation is performed on a 
256-site lattice with 50 random realizations. The small devi¬ 
ation in the end is a finite size effect when the subsystem size 
L approaching to the half-system size, i.e. log 2 128 = 7. 


Fig.[^ shows our numerical result for the quantum 
Ising model in Eq. 0 with coupling coefficients drawn 
from uniform initial distributions (at Tq = 1). At the 
marginal MBL criticality (A, B and E), the EE follows 
the logarithmic scaling: 


SsiL) = 3 InA. 


(29) 


where c' is the effective central charge at the strong dis¬ 
order fixed-point, which is related to the central ch arge 
c of the same system in the clean limit by d = cln2.ff2£l 
The point A(l, 0,1) corresponds to a free Majorana chain 
with c = 1/2. The point 5(1,1,0) corresponds to two 
free Majorana chains with c = 1. The point 5(2,1,1) 
is in between A and 5, and corresponds to an inter¬ 
acting Majorana chain with c = 1/2, because hi is a 
relevant perturbation which always drives the two Majo¬ 
rana chains into one. The corresponding effective central 
charges d in all these cases agree with our numerical re¬ 
sults as shown in Fig.fTHl 

In the MBL phase (C, F and G), the EE saturates 
to the area law (which is a constant for the ID system). 
The point G(l, 0,0) is deep in the SG phase, so its EE is 
exactly Ibit {Se = In 2) as expected for the Greenberger- 
Horne-Zeilinger (GHZ)i^ state of the large block spin (or 
from the Majorana zero modes at the entanglement cuts). 
The point G(0,1,1) is deep in the PM phase, so its EE 
must vanish (S”^ = 0) because all the eigenstates are on¬ 
site (atomic) direct product states. The point 5(1,1,1) is 
a generic point of the interacting Majorana chain, whose 
EE will initially grow a little with L and quickly saturate 
to the area law at a non-universal value. We can see that 
all these cases can be handled by SBRG nicely. 





FIG. 17: Phase diagram of the quantum Ising model Eq. (111 
in ternary graph of (Jo, Ao, ho) = {Jo, Ko, ho)^^^° ■ The color 
plot shows the half-system-size EE Se ■ The black curve traces 
out the SG-PM phase boundary following the local maxi¬ 
mal of Se- The colored points are labeled by their coordi¬ 
nates (Jo, Ao, ho) as A(l, 0,1), B(l, 1, 0), G(l, 0, 0), D(4,1,1), 
5(2,1,1), 5(1,1,1), G(0,1,1). The calculation is performed 
on a 256-site lattice. 


Because the EE diverges logarithmically along the 
phase boundary between the two MBL phases, it can 
be used to map out the phase diagram. We calculate 
the half-system-size EE on a sufficiently large system 
(256-site), and the result is shown in Fig. 17 The phase 
boundary is traced out along the local maximum of the 
EE in the phase diagram. It was also proposed in Ref.lMI 
to detect the phase boundary by the divergent EE fluc¬ 
tuation. We also tried that method as well and obtained 
basically the same phase diagram (not shown). 
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D. Locality of the Effective Hamiltonian 


At the end of SBRG flow, we arrived at the effective 
Hamiltonian of the emergent qubits in the holo¬ 
graphic bulk. 

HeS = '^eiTi-[-'^tijTiTj+ ^ eijkTiTjTk-\ -, (30) 

i i<j i<j<k 


which can be considered as the RG transform of the ini¬ 
tial physical Hamiltonian H to the 2nd order precision, 
as Hes — U^qHU^g- So how does i?eff look like in the 
holographic bulk? Is is still a local Hamiltonian? To an¬ 
swer these questions, we must first assign the position to 
each emergent qubit in the holographic bulk. 

The position of each emergent qubit Ti can be labeled 
by (xi,€i), where Xi is the real-spa ce coo rdinate of the 
qubit in the Glifford circuit in Fig.|9pi^ and is the 
energy scale associated to Ti which is taken from the 
single-body terms in TJeff- It is reasonable to use the 
qubit position in the Glifford circuit directly as the real- 
space coordinate, because the Glifford circuit preserves 
the real-space locality from the holographic boundary, as 
it only contains at most ~ InL layers of local Clifford 
gates. If one wishes to be more precise about the real- 
space coordinate, one can map the emergent qubit back 
to the stabilizer Ti on the holographic boundary, and use 
certain center-of-mass coordinate of the stabilizer Ti as 
the real-space coordinate for the emergent qubit. How¬ 
ever center-of-mass method gives roughly the same co¬ 
ordinate, because the qubit position method already en¬ 
sures Xi to fall in the real-space support of fi and we 
have demonstrated that fi has nice locality in the MBL 
system (in Fig. 14). 


uj" 1 (a) . . * 

— ^ • t* • 




FIG. 18: Effective Hamiltonian Heff in the holographic bulk 
for (a) {Jo,Ko,ho) = (1,1,1) and (b) {Jo,Ko,ho) = (2,1,1) 
on a 32-site lattice with periodic boundary condition. Each 
black dot denotes an emergent qubit Ti placed in the holo¬ 
graphic bulk according to its energy scale Si. The 2-body 
(3-body) terms in the effective Hamiltonian are represented 
by the links (triangles) between the emergent qubits. The 
darker object corresponds to the stronger interaction term in 
Heft. 


As we have pinned down the positions of the emergent 


qubits Ti, we can draw pictures of the effective Hamilto¬ 
nian ifeff in the holographic bulk. Two typical examples 
are shown in Fig. |18[ where the 2-body and 3-body terms 
are shown respectively as the links and triangles among 
Ti- Fig.[I^a) shows an example in the MBL phase, which 
demonstrates nice real-space locality: the links and tri¬ 
angles are typically small, and there is no long object 
stretching through out the system. Fig.[T8|(b) is an ex¬ 
ample at the marginal MBL criticality, where the en¬ 
ergy scale is extended much deeper into the IR limit, and 
longer links/triangles are also observed, which is consis¬ 
tent with the delocalization at the criticality. In both 
pictures, we can see the strong links and strong triangles 
tend to appear together in the same region, which im¬ 
plies that there should be some correlation of the 2-body 
terms with the 3-body and higher many-body terms in 
ifeff, as pointed out in Ref. [80] in a different context. 




FIG. 19: Norm of the 2-body energy coefficients eij as a func¬ 
tion of the distance Xij between the emergent qubits for (a) 
(Jo, Ko, ho) = (1,1,1) in log-linear plot and (b) (Jo, Kq, ho) = 
(2,1,1) in log-log plot. The calculation is performed on a 32- 
site lattice with 1000 random realizations. 


Now we will focus on the 2-body energy coefficients etj, 
and study their locality quantitatively. We classify the 
tijTiTj terms in the effective Hamiltonian iJeff by their 
real-space distances \xi — Xj\, and calculate their Frobe- 
nius norm as a function of the distance. 



It is verified that the norm decays exponentially 

(32) 

with the real-space distance Xij in the MBL phase as 
shown in Fig.[l^a), which was also proposed in Ref.l^Tl- 
[2H1 The real-space decay of ||eijj| is related to the loga¬ 
rithmic growth of EE after a global quench in the MBL 
systemas S{t) = Soo^ln(Jot), which is an MBL in¬ 
trinsic phenomenon that is not present in the Anderson 
localization. At the marginal MBL criticality, the norm 
follows a power-law behavior ||eijH ^ as is shown in 
Fig.[^b). 
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E. Holographic Hamiltonian 


To further study the entanglement structure among 
the emergent qubits, we map the initial physical Hamil¬ 
tonian from the holographic boundary into the bulk by 
the Clifford circuit, and define the bulk Hamiltonian as 
the holographic Hamiltonian -Hhoi, 

ilhoi = (33) 

The holographic Hamiltonian iJhoi acts on the emergent 
qubits in the holographic bulk. Because the unitary 
transform Uc\ can be carried out exactly, i7hoi has ex¬ 
actly the same spectrum as H, and its energy eigenstates 
are the corresponding exact eigenstates of H mapped into 
the holographic bulk by the Clifford circuit Uqi- 

Because Uci does not exactly diagonalize the Hamilto¬ 
nian i7, i7hoi contains off-diagonal terms, which give rise 
to the entanglement among the emergent qubits in the 
bulk eigenstates of Hhoi- The portion of the off-diagonal 
terms in the holographic Hamiltonian Hhoi can be mea¬ 
sured from the following ratio: 


Tr(ghoi - diagHhoi)^ 


= S:. (34) 


where diagiJhoi = E{ri} l{'^i})({'^JI^hoi|{'ri})({Ti}| de¬ 
notes diagonal part of i7hoi, and on the right-hand- 
side, 5E‘^ and E"^ are respectively the same as E q.([^ 
and Eq. (25). We have checked this ratio 5E'^/E'^ in 
Eig.[TTJ and verified that the off-diagonal term does not 
dominate the holographic Hamiltonian i7hoi for MBL sys¬ 
tems. 

To discover the geometry in the holographic bulk, one 
can in principle diagonalize i7hoi and study the mutual 
information lij between two emergent qubits Ti and tj 
on the eigenstate, 


+ (35) 

where Si (Sj) is the EE of each single qubit, and Sij is the 
EE for both qubits. Then the distance dij between two 
qubits can be defined following the proposal in Ref. l81l82l 

d^j = (36) 

to 

where Iq = 2 In 2 is the maximal mutual information 
between two qubits. In the strong disorder limit, the 
holographic Hamiltonian iJhoi becomes diagonal, and its 
diagonal terms also coincide with the effective Hamilto¬ 
nian HeS- In this limit, the emergent qubits are disen¬ 
tangled from each other, so they are infinitely far from 
each other, and the bulk space is fragmented into isolated 
points.l^ Away form the strong disorder limit, the off- 
diagonal terms in ifhoi start to grow (see Fig. [n]) , which 
leads to the resonance between the emergent conserved 
quantities and the entanglement among the bulk qubits. 
From the geometry perspective, the emergent qubits are 
getting closer to each other, as if the space is contracting. 



FIG. 20: Pictures of the off-diagonal terms in the holo¬ 
graphic Hamiltonian Hhoi for (a) {Jo, Kq, ho) = (1,1,1) and 
(b) {Jo, Kq, ho) = (2,1,1) on a 32-site lattice with periodic 
boundary condition. Each black dot denotes an emergent 
qubit Ti in the holographic bulk. The off-diagonal terms of 
2-qubit (3-qubit) resonances are represented by the links (tri¬ 
angles) between the emergent qubits. The darker object cor¬ 
responds to the stronger off-diagonal term in Hhoi- 


To gain more intuition, we draw the off-diagonal terms 
of Hhoi in the holographic bulk for both the MBL 
Fig.pU|(a) and the marginal MBL Fig.[20|(b) systems. The 
off-diagonal terms create the spacial connectivity in the 
holographic bulk and drives the emergent qubits closer 
to each other. In Fig.[20](b) one can also observe a darker 
resonant clusteil^ around Xi ~ 16, in which the qubits 
are connected by more and stronger off-diagonal terms. 
We suspect that such regions are closer to thermalization. 
It is conceivable that if iLhoi is dominated by the off- 
diagonal terms in a local region of the holographic bulk, 
the Hamiltonian in that region will look like a random 
matrix, and the corresponding eigenstates will look like 
random states,!^ therefore the region is locally thermal- 
ized. In the thermalized region, the geometry becomes 
non-local as every qubit is almost maximally entangled 
with the rest of the qubits in the region. Therefore the 
thermalized region can also b e viewe d as a small black- 
hole in the holographic bulk.ESHini! jn the MBL-ETH 
transition, small blackholes will first emerge from the UV 
region, and then merge into larger blackholes. When a 
large blackhole covers the IR region, the whole system 
will be thermalized, since the low- freq uency physics is 
then dominated by the ETH phase.^^ The above bulk 
scenario is dual to ergodic puddle percolation on the holo¬ 
graphic boundary, which was proposed in Ref. l34l89l90l 
SBRG approach will break down in the local ther¬ 
malized region, where the holographic Hamiltonian Hhoi 
is dominated by strong off-diagonal resonances among 
the emergent qubits. One key assumption of SBRG is 
that there exists a set of local quasi-conserved quantities 
in the MBL system which will not resonate with each 
other due to their different energy scales, so that the 
off-diagonal terms can be treated by non-degenerate per¬ 
turbation. Thus if we found that the off-diagonal terms 
actually dominate Hhoi after SBRG calculation, it would 
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imply that the assumption fails and the RG result is 
not reliable. In that case, the emergent qubits found 
by SBRG in the local thermalized region are no longer 
quasi-conserved quantities. 

Towards thermalization, SBRG becomes not only in¬ 
accurate, but also inefficient. Consider a local thermal¬ 
ized region, in which the Hamiltonian H is dominated 
by the off-diagonal terms, the 2nd order perturbation 
Hq + 'S —>■ Hq — in Eq. ([^ will roughly square 

the number of terms in H at each RG step. As a result, 
the number of Hamiltonian terms in the local thermal¬ 
ized region will grow super-exponentially with the RG 
flow. The computation complexity of each SBRG step is 
proportional to the number of Hamiltonian terms, so the 
RG will quickly get stuck at the boundary of the local 
thermalized region (also the horizon of the small black- 
hole) due to the fast growth of the complexity. Although 
SBRG can not approach the ETH phase, we can still 
observe indications of local thermalizations in the MBL 
phase from the increasing off-diagonal resonances in TJhoi 
or from the growth of the Hamiltonian terms along the 
RG flow. 


V. SUMMARY AND DISCUSSION 

In summary, we have demonstrated that SBRG is 
an efficient numerical method to study the MBL (or 
marginal MBL) systems. The full energy spectrum and 
all the many-body eigenstates of an MBL system can be 
approximately found with an algorithm complexity that 
scales linearly with the system size, which is much more 
efficient than the ED method which has an exponential 
complexity. Although we have only explored ID mod¬ 
els in this work as examples, the general formalism of 
SBRG is readily applicable to higher dimensional mod¬ 
els on any lattice. The accuracy of SBRG is controlled 
by the randomness of the system. We have shown that, 
deep in the MBL phase, SBRG flows to the strong disor¬ 
der limit and becomes asymptotically exact. However, as 
we tune the system towards the ETH phase, SBRG will 
be more and more inaccurate and inefficient. Although 
SBRG fails to approach the ETH phase and the ther¬ 
malization transition, we can still use it to identify small 
locally thermalized regions in the MBL phase. Neverthe¬ 
less, the MBL-ETH transition can be described by other 
RG sche mes be yond SBRG, such as the phenomenologi¬ 
cal RSRGP^M! based on merging the locally thermalized 
regions. Ref.EH further points out that both the MBL 
and the ETH states can be modeled on the classical level 
by Clifford circuits. It is interesting to note that the 
Clifford circuit is naturally generated by SBRG as the 
classical approximation of the EHM circuit to describe 
the MBL states. Whether it is possible to design the RG 
scheme that generates the ergodic Clifford circuit for the 
ETH state is another interesting problem to explore. 

As a Hilbert-space-preserving RG scheme, SBRG can 
also be interpreted as a realization of the holographic 


duality. More precisely, an Clifford circuit can be con¬ 
structed from SBRG flow, which maps the MBL system 
from the holographic boundary to the bulk. It turns 
out that the holographic bulk degrees of freedom are 
the emergent conserved quantities of the MBL system 
which are (approximately) governed by the MBL fixed- 

point Hamiltonian iJeS = + J2ij H-■ In 

the strong disorder limit, the degrees of freedoms in the 
holographic bulk are disentangled, and the spacial geom¬ 
etry is fragmented. Away from the strong disorder limit, 
we can show that the locally-thermalized regions emerge 
in the MBL state, which may be interpreted as the small 
blackhole formation in the holographic bulk. However, 
under the present SBRG scheme, the emergent locality 
along the energy scale can not be observed, and very 
strong UV-IR mixing prevails in the holographic bulk. 
We think it is because that each bulk qubit is renormal¬ 
ized independently and abruptly in one RG step, such 
that the entanglement between the UV and the IR qubits 
can not be efficiently removed. It will be desirable to im¬ 
prove the RG scheme such that all qubits are renormal¬ 
ized jointly and the Hamiltonian is diagonalized gradu¬ 
ally by small unitary transforms. 

Note: During the completion of t his manuscript, we 
become aware of a similar worli^dJ^ by L. Rademaker, 
which introduces a different RG scheme that also finds 
the emergent conserved quantities of the MBL systems by 
consecutive unitary transforms of the Hamiltonian. The 
difference is that we keep all orders of interaction terms 
but treat the off-diagonal resonance perturbatively, while 
there the resonance is solved exactly but the higher order 
interactions are truncated. 
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Appendix A: Technical Details of SBRG 
1. Clifford Group Rotations 


• If [/x| contains the index 1 or 2, such as or 

(i^ = 0,1, 2,3): rotate that index to 3 while 
eliminating other index to 0 by a (74 rotation as 


For our purpose, we generate the Clifford group by the 
(74 rotation i?c 4 (like f phase gate) and the swap gate 
-Rswap- The (74 rotation is defined by its generator 
as 


= exp 


Its adjoint action on a matrix is given by 

_ f if commute, 

“ I icrt'^lcrl'^1 if crM, crl'^l anti-commute. 


(A2) 


The swap gate is specified by i O j, where i and j label 
the two qubits to be exchanged. For example, on the 
two-qubit level 

7?swap(1 O 2) = (A3) 


which can be generalized to any pair of qubits straightfor¬ 
wardly. The adjoint action of the swap gate i?swAp(* ^ 
j) on a matrix is simply exchanging the indices 
O /ij in [/i]. One can see that both the rotation 
and the swap gate can be implemented efficiently on the 
algebraic level by manipulating the Pauli indices of the 
matrices. In SBRG algorithm, the Hamiltonian is never 
spelt out as a matrix (not even as a sparse matrix), all we 
kept is a table of coefficients (the coefficient in front 
of the matrix tjl^^). The Clifford group rotations is im¬ 
plemented on the index [/i] (with a possible sign change 
of /i[^] when needed)^ _ 

According to Eq. ( |A2[ ) , if we want to bring a matrix 
to another matrix crt^l which anti-commutes with 
the original one, we just need to perform one Ca rotation 

i?C4(iCT[^lcrH): 


a 


p.[A] 


(A4) 


However, if we want to bring a matrix to a com¬ 
muting matrix we can find an intermediate matrix 
cr[“l which anti-commutes with both and and 
perform two consecutive C 4 rotations, which we called a 
double-(74 rotation. 


i?C4(i' 


rH„[- 


-)► cr^ 






->■ a 


[a«] 


(A5) 


Finally if the certain Pauli indices of a matrix do not 
appear at the position that we want, we can use the swap 
gate to permute the indices to the intended position. 

With these, we come up with the following protocol 
to rotate (A = 0, 3 and /r = 0,1, 2,3) to the block 

diagonal form 




(A6) 


then move this index 3 to the intent qubit position 
by a swap gate if necessary. 

• Else if [^] has no index of 1 or 2, but the index 3 
exists: 


— If there are multiple indices of 3 in [A] [/i]: we 
need to perform the double-(74 rotation: 

* If there is an index 3 right at the intent 
qubit position, such as {k = 0,3): 

the double-C 4 rotation goes as 


^[A]3[k] ^[o...]l[o...] 

flc4(-<T[°-l"t°- l)^ ^[0-]3[0-], 


(A7) 


* Else the intent qubit position must host 
the index 0: then find the last qubit po¬ 
sition of index 3, and perform the above 
double-(74 rotation at that position, and 
finally use a swap gate to bring that qubit 
to the intent position. 

— Else there is only one index 3 in [/i]: just use 
a swap gate to bring that qubit to the intent 
position. 

• Else [/r] is [0 • • • ]: no Clifford group rotation needed. 
In fact, such matrix crl^ko -] gj^jould have already 
been ascribed to the effective Hamiltonian and 
should not actually appear in the residual Hamil¬ 
tonian. 

At the end of SBRG flow, the Hamiltonian is diag¬ 
onalized, and all the emergent conserved quantities are 
identihed as the emergent qubits in the holographic bulk. 
Each emergent qubit Ti is positioned by both its real- 
space coordinate Xi and its energy scale e^. Following the 
above protocol, the emergent qubits Ti will be arranged 
(roughly) by their energy scales e^, i.e. ei ^ £2 ^ ■ 

This has the advantage of improving the algorithm ef¬ 
ficiency, because the physical and the emergent Hilbert 
spaces can be easily separated. However the real-space 
locality is lost, as the real-space coordinates Xi are scram¬ 
bled by the swap gates. For many physical applications, 
it is more desired to reorder the emergent qubits by 
their real-space coordinates Xi. To restore the real-space 
locality, we simply need to push all the swap gates out 
of the Clifford circuit, and apply them to the emergent 
qubits to bring the qubits back to the real-space ordering. 
When the swap gate passes through a C 4 gate, the Pauli 
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indices of the C 4 generator should be swapped accord¬ 
ingly. In this paper, we always assume that the real-space 
locality is restored after SBRG flow, and are ordered 
by their real-space coordinates Xi = i. Under this imple¬ 
mentation, the Clifford circuit will only contain the C 4 
rotations and no swap gate. 


2. SchriefFer-Wolff Transformations 


Consider the Hamiltonian H = iJo -I- A -|- E, where 
i?o = is the leading energy scale in the block- 

diagonal form, A is the rest of the diagonal block 

E (A 8 ) 

A=0,3 [^t] 

and E is the off-diagonal block 

>^=1,2 M 

The matrices i?o, A and E are all Hermitian, and satisfy 
Hi = hi, HqA = AHo, i?oS = -E7?o. (AlO) 


We can eliminate the off-diagonal block to the 2nd order 
in by the SchriefFer-Wolff transformation H —>■ 

S'^HS, where 




^ 2hl^°^ 8hl^^ 


(All) 


Using the properties in Eq. (AlO), it is straightforward 
to verify that to the 2 nd order perturbation we have 


H S'^HS 

1 T 1 , , (A12) 

~ iJo + A-I-^^i7oE +-^Ho[E,A\. 

The first three terms are block-diagonal, and the last 
term is block-off-diagonal, which can be checked by ex- 
aming their commutation relations with Hq. We simply 
project out the block-off-diagonal part i/o[S, A]/(2h|). 
Although it is a 2nd order term, but if we further ro¬ 
tate it to the diagonal block, it will only produce a 4th 
order correction to the diagonal block, which is negligi¬ 
ble. Therefore we end up with the 2nd order effective 
Hamiltonian in Eq. ([^ after the SchriefFer-Wolff trans¬ 
formation. 

Since the 2nd order perturbation iJoE^/(2/i§) will gen¬ 
erate many terms and cause the Hamiltonian to grow. 
So in order to prevent the uncontrolled growth of the 
Hamiltonian that may crash SBRC program on the com¬ 
puter, we need to control the growth rate by truncating 
the small terms generated in the 2 nd order perturbation. 
In practice, we first set a maximum growth rate r (say 


r = 2). In each RC step, we count the number of terms in 
E and denote it as A/e, then we only keep the leading cA/e 
terms in HqT,^/{ 2 hl) and discard the rest of the smaller 
terms. In this way, the Hamiltonian will not grow too 
fast, and the complexity of the algorithm is controlled. 
The truncated terms can be collected for the error esti¬ 
mation afterward. We found that deep in the MBL phase 
very few terms are truncated, but around the marginal 
MBL critically more terms will be truncated. One can 
also adjust the maximum growth rate r to see if SBRC 
result converges in the r —>■ 00 limit. We found that 
in most of the cases, r = 2 is already good enough to 
converge the Clifford circuit generated by SBRC. 


Appendix B: Entanglement Entropy of Stabilizer 
States 


Let 6 = {fi} be the set of emergent conserved quanti¬ 
ties represented in the original physical Hilbert space (as 
Pauli operators). Then every energy eigenstate in the 
MBL spectrum can be approximated by a state 
stabilized by 6 as Vz : TzI'I'it-.}) = In the 

quantum information terminology, © is a complete set of 
stabilizers, and |'L{ri}) is a stabilizer state. 

The bipartite EE of a stabilizer state can be calcu¬ 
lated by the highly efficient method developed in Hef. 11051 
which will be briefly reviewed here. Suppose we are in¬ 
terested in the EE Sa in a subsystem A of the state 


'S'a =-TrpAlogsPA, 

PA = Tr^ |^'{n})(^{n}l, 

where A denotes the complement of A. We can first clas¬ 
sify the stabilizers fi into three sets: the stabilizers only 
supported in A as ©.4 = Tr^ ©, the stabilizers only sup¬ 
ported in A as &A = Tr^ ©, and the rest of the stabilizers 
supported in both A and A as &aA = © — ©A — ©A- 

The EE is fully determined by &aa- Because ©aa 
contains all the stabilizers that would be broken by the 
entanglement cut, and thus can not be used to stabilize 
the state in A (or in A), which gives rise to the entropy 
as the state in the subsystem can not be fully deter¬ 
mined. Naively the EE would just be proportional to 
the number of stabilizers in ©aA) however this idea must 
be refined, because some stabilizers in ©aa actually 
hidden stabilizers which can be localized to either A or 
A by multiplication with other stabilizers, and hence not 
contributing to the EE. 

To reveal the hidden stabilizers in ©aai can rewrite 
each stabilizer fi £ ©aa as a direct product ii = 
explicitly, where is the part of the Pauli operator 

Ti that is supported only in A (A). Then rf may not be 
commuting with each other, and the statement is that 
the hidden stabilizers are in one-to-one correspondence 
to the center of the Pauli group generated by ff-. 
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To reveal the center, we can first construct the matrix 
M of anti-commutativity among the operators 


yTi,fj € &aa ■ 


Mij — 


1 = -TfT^, 

0 otherwise. 


(B2) 


will not contribute to the EE, so the nullity must be 
subtracted, and the EE (in unit of bit) is given by the 
Z 2 -rank of M 


Sa 


1 

2 


rankM (over Z 2 ). 


(B3) 


It can be shown that the center corresponds to the null 
space (the kernel) of the integer matrix M modulo 2. 
Because the null space basis are hidden stabilizers that 


The Z 2 -rank of a matrix can be found using Gaussian 
elimination, and the EE can be simply calculated from 
the algebraic structure of the stabilizers. 
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Originally EHM is short for “exact” holographic mapping. 
However, in this paper, EHM circuit is actually an approx¬ 
imation of the exact diagonalization circuit, so to avoid 
confusion of the mean of “exact”, we change the name to 
“entanglement” holographic mapping. 

In higher spacial dimension, the fermion model will be 
mapped to the non-local qubit model in general. Although 
SBRG algorithm complexity may grow exponentially in 
the absence of locality, we may still consider the non¬ 
local qubit model as our most generic stating point for 
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theoretical discussions. 

The coefficient here is slightly different from the 
used previously, in that the 2nd order correction to the 
- ] been absorbed into the new coefficient 

here. 


For comparison, 5E^ /E'^ on random states will be close 
to 1. 

In general Xi can be a real-space coordinate vector for 
higher dimensional problems. 



